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ABSTRACT 

A two dimensional steady-state Navier-Stokes 
solver has been upgraded to include the effects of 
frozen and equilibrium air chemistry for applications 
to high speed flight vehicles. To provide a computa- 
tionally economical first order approximation to the 
high temperature physics, variable thermodynamic 
data is used for the chemically frozen mode to allow 
for a variation with temperature of the air specific 
heats and enthalpy. For calculations involving air in 
chemical equilibrium, a specially modified version of 
the NASA Lewis Chemical Equilibrium Code, CEC, is 
used to compute the chemical composition and resul- 
tant thermo-chemical properties. The upgraded solver 
is demonstrated by comparing results from calorically 
perfect (C p =constant), thermally perfect (frozen) and 
equilibrium air calculations for a variety of geometries, 
and flight Mach numbers. 

INTRODUCTION 

The two dimensional steady state Navier-Stokes 
solver that has been upgraded to account for first 
order real gas effects is the PARC2D code. The 
PARC2D code is a highly modified version of AIR2D 1 
that has been made into a practical engineering anal- 
ysis tool for propulsion oriented problems by Cooper, 
et. al. 2 of Sverdrup Technology at AEDC. PARC 
solves the Reynolds averaged Navier-Stokes equa- 
tions using the approximate factorization algorithm 
of Beam and Warming 3 . To avoid the expense of 
inverting the large block banded matrix system of 
AIR2D, Pulliam 4 diagonalized the implicit operators, 
which resulted in signifigant improvements in compu- 
tation time. Artificial dissipation was also added as 
proposed by Jameson 6 to improve convergence and 
shock capturing capability, which then resulted in the 
ARC2D series of codes. This code was then modi- 
fied by Cooper et. al, 1 for use in solving propulsion 
oriented problems at the Arnold Engineering Devel- 
opment Center. These modifications resulted in the 
PARC series of codes. 


The PARC codes are well validated and accepted 
and are in wide use by the aerodynamics community. 
Through the use of a unique method of grid patch- 
ing and the segmented application of boundary con- 
ditions, the PARC series of codes are very useful for 
solving propulsion oriented problems. This solver was 
originally developed to calculate flows of a single con- 
stituent that is assumed to be a calorically perfect 
gas. This assumption is usually good for lower speed 
flows, but for high speed calculations with resultant 
high enthalpies, wall temperatures and heat transfer 
rates can be overpredicted. To permit more realistic 
calculations of high speed flows, this useful flow solver 
has been modified to account for flows of thermally 
perfect gases and flows where air is in thermal and 
chemical equilibrium. 

ANALYSIS 

PARC solves the Reynolds averaged Navier- 
Stokes equations on a curvilinear coordinate system 
written in divergence form, given below 

dtq + d x E + d v F = Re~ 1 (d x R + d y S) (l) 

where q is a vector containing the conservation vari- 
ables: 
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The vectors E and F are the inviscid flux vectors: 
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where r xx = (A + 2/i)u x + Xv v , r xy = n[u y + v x ), and 
T y y = (A + 2 fl)v y + A U*. 


The original formulation assumed the gas to be 
calorically perfect and thermally perfect 2 . This as- 
sumption results in a simple algebraic relationship be- 
tween the pressure, temperature and flux variables, 
namely 

( 6 ) 

r = l( 7 -,)(?i_iM^l>) ( 7 ) 

?i 1 

where 7 is constant. For a gas where the specific heats 
are variable and the gas follows the ideal equation of 
state, the pressure and temperature are interrelated 
through the definition of the enthalpy, 


h(T) = E yi hi{T) 
f T 

= / C p dT+Ah° f 
J Tr 

= e+ £ ( 8 ) 

where y; is the mass fraction of the i-th constituent of 
the gas. By assuming the ideal gas equation of state, 
eqn. (8) is manipulated and written in terms of the 
flux variables as 


Mr)-Jir=2i-l[(a) ! + (ffi) ! ] (9) 

For this analysis, the thermodynamic data for C p and 
h are in the form of polynomials in temperature for 
each species of the ideal gas mixture. These data have 
been calculated by statistical mechanics for the tem- 
perature range of 200 to 15000 Kelvin for the species 
N, NO, N2, O, N, 02, e-, N+, NO+, O+, N2+ and 
02+ 6,7 ' 8,9,ic> . The functional form for the specific 
heat and enthalpy of each species is 

. = + 02^ + 03T 2 + a^T 3 + asT* (10) 
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The polynomials use two sets of coefficients valid on 
the ranges of 200>T>5000 and 5000>T>15000 Kelvin 
and are constrained to match at the midpoint temper- 
ature of 5000 Kelvin. Once the chemical composition 
is known, eqn. (9) is solved iteratively for the temper- 
ature using a Newton-Raphson iteration, which then 
yields the pressure through the ideal gas equation of 
state. To calculate the derivatives of pressure needed 
in the flux jacobians an effective gamma approach is 
used similar to the methods presented in references 
11 and 12. The boundary condtion routines for slip, 
no-slip adiabatic and specified temperature walls were 
modified to be consistent with the thermodynamics 
and resulting equation of state. For frozen flows, the 
chemical composition is held constant throughout the 
flowfield at the standard air composition, while for 
the equilibrium air calculations the species solutions 
are obtained separately using an equilibrium solver. 

CHEMICAL EQUILIBRIUM SOLVER 

, The species and thermodynamic solutions needed 
for the equilibrium air calculations are obtained in 
a "black box” fashion using the newly developed 
equilibrium air package developed at NASA Lewis 
by McBride et. al. 6 . This package is based on a 
Helmholtz free energy minimization procedure fully 
described by Gordon and McBride in reference 13. 
This package requires only the local density and inter- 
nal energy (which are simply related to the flux vari- 
ables) as input, and returns as output the resultant 
pressure, temperature, equilibrium speed of sound 7, 
and chemical composition. The chemical composition 
is then used in a consistent manner to compute the 
needed thermodynamic properties. A special proce- 
dure is used to initialize the chemical equilibrium cal- 
culations as well as a procedure to determine freezing 
conditions. The development of this unique equilib- 
rium solver is still underway. A more complete de- 
scription of the package, procedures and governing re- 
lations is available in reference 6. 

RESULTS 

To demonstrate the new capability of this solver, 
calculations of a variety of geometries and flight con- 
ditions are presented. Comparisons are made between 
the calorically perfect, thermally perfect and chemical 
equilibrium air calculations for these hypersonic flows. 
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(1) Moo = 7 Hypersonic Ramp Case 

Mach 7 flow over a 15° ramp is calculated at a 
freestream pressure of 1 atmosphere and temperature 
of 350 Kelvin. To highlight the real gas upgrades an 
adiabiatic wall boundary condition was chosen for this 
laminar calculation. A simple sheared grid was used 
with clustering along the lower wall for this 100 (axi- 
ally) by 50 (normally) grid. The upper boundary was 
treated as a symmetry plane. Figure 1 shows the com- 
puted pressure contours for the calorically perfect cal- 
culation while Figures 2 and 3 compare the wall pres- 
sures and temperatures (respectively) for the three dif- 
ferent thermo-chemistry modes. Figure 4 indicates the 
amount of dissociation calculated in the equilibrium 
calculation by showing the amount of NO along the 
wall. The L2 norm of the residuals was driven down 
by over 4 orders of magnitude in approximately 2000 
iterations for the three different calculations. 

The pressure contours in Figure 1 show the lead- 
ing edge shock merging with the ramp shocks, and the 
shock /shock interraction of the upper and lower ramp 
shocks. Examination of Figures 2 and 3 shows marked 
differences in both wall pressure and temperature be- 
tween the non-real gas and the real gas calculations. 
The calorically perfect calculation (C p =constant) in- 
dicates a small separation/interraction region extend- 
ing approximately 0.2 meters ahead of the ramp turn 
(the ramp turn is located at X=0.3 meters). This sep- 
aration/interraction region is shortened considerably 
in the frozen air (variable C p ) and equilibrium air cal- 
culations. It appears that by accounting for more ther- 
modynamic degrees of freedom a larger turning angle 
could be accepted by the flow before separating. Af- 
ter the initial temperature and pressure rise caused by 
the ramp, it is interesting to notice that the frozen air 
and equilibrium air wall temperatures are fairly close, 
while the calorically perfect wall temperature is much 
higher for this lower Mach number case. 

(2) Moo = 10 Hypersonic Ramp Case 

A Mach 10, 10° ramp case was calculated with 
freestream conditions of 1 atmosphere pressure and 
300 Kelvin temperature. This case was also run lam- 
inar and with an adiabatic wall to highlight the real 
gas effects. The same grid dimensions and near wall 
stretching used in the Mach 7 case were used for this 
calculation. The L2 norm of the residuals was driven 
down by over 4 orders of magnitude in approximately 
4000 iterations for all three different thermo-chemistry 
modes. 

Figure 5 illustrates the computed pressure con- 
tours for the calorically perfect calculation, while Fig- 
ures 6 and 7 show the wall pressure and tempera- 
ture comparison between the three thermo-chemistry 
modes. The extent of air dissociation is indicated in 


Figure 8 which shows the computed mass fractions of 
NO along the wall. Examination of the pressure and 
temperature plots indicates a fairly large separation 
region ahead of the ramp for the calorically perfect 
calculation with a milder pressure rise profile. The 
frozen and equilibrium calculations show a progres- 
sively smaller interraction region with a sharper pres- 
sure rise. In this case, the final wall temperatures 
indicate a much larger difference between the calor- 
ically perfect, frozen air and equilibrium air calcula- 
tions. The amount of dissociation is nearly twice that 
of the Mach 7 case which is indicated by the fairly 
high amount of NO present at the wall. 

(3) P2 Hypersonic Inlet at Moo = 10. H=100 Kft. 

The P2 inlet geometry, taken from reference 14, 
was used to represent a large scale inlet at a flight 
Mach number of 10 at an altitude of 100,000 feet. For 
this Mach number and altitude the Reynolds num- 
ber is low enough to assume that the flow is initially 
laminar, and to remove any issues about transition 
to turbulence and the correct modelling of the turbu- 
lence, the flow is considered to remain laminar along 
the length of the inlet. It is also assumed that the wall 
temperature remains constant at 1000 Kelvin. A 100 
by 80 grid was used with near wall clustering on both 
the ramp and cowl walls. Figure 9 shows the pres- 
sure contours calculated in the calorically perfect gas 
mode. The comparisons between the different thermo- 
chemistry calcultions for the ramp and cowl pressures 
and temperatures are shown in Figures 10 and 11, 
respectively. Figure 12 compares the computed exit 
plane temperature profiles between the three different 
calculations. For the calorically perfect and frozen air 
calculations the L2 norm of the residuals was driven 
down by over^ orders of magnitude in approximately 
5000 iterations while the equilibrium calculation resid- 
ual was driven down only 3 orders. 

Examination of the ramp and cowl pressure com- 
parisons (Figures 10 and 11) shows little difference be- 
tween the different thermo-chemistry modes. The exit 
plane temperature profiles show minor differences in 
the temperature peaks and in the core flow regions. 
These results show that for this relatively low hyper- 
sonic Mach number and for the "cooled” wall condi- 
tions that the real gas effects for this case are negligi- 
ble. 

SUMMARY AND CONCLUSIONS 

A two dimensional steady state Navier) Stokes 
flow solver has been modified to account for real gas 
effects and has been demonstrated on various geome- 
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tries at hypersonic Mach numbers. The chemically 
frozen air modifications account for the variation of 
the air thermodynamics using polynomials in temper- 
ature for the specific heats and enthalpies that are 
valid over the temperature range of 200 to 15,000 
Kelvin. A specialized version of the NASA Lewis 
Chemical Equilibrium Code, CEC, was used to com- 
pute the chemical composition and resulting thermo- 
dynamic properties given a thermodynamic state spec- 
ified by the flow code flux variables. 

Three separate geometries and free stream Mach 
numbers and conditions were calculated and com- 
parisons were made between the calorically perfect 
(C p =constant), frozen air (variable C p ) and equilib- 
rium air thermo-chemistry modes. Considerable dif- 
ferences were shown in wall pressures and tempera- 
tures for the two adiabatic wall hypersonic ramp cases. 
For these two cases, the adiabatic wall boundary con- 
ditions enabled the flow near the wall to reach near 
stagnation conditions, resulting in signifigant differ- 
ences in the computed wall temperatures due to dis- 
sociation in the equilibrium calculations and to the 
variatiqn of C p in the frozen air calculations. The 
large scale inlet simulation showed negligible differ- 
ence between the different thermo-chemistry modes. 
This is believed to be due to the specification of wall 
temperatures at levels low enough not to cause suffi- 
cient dissociation of the air. Although it cannot be 
extrapolated to higher Mach numbers or lower alti- 
tudes (resulting in higher stagnation conditions), it 
appears that for the large scale inlet case computed 
in this study little difference can be noted between 
the calorically perfect, frozen air and equilibrium air 
results. 

ACKNOWLEDGEMENTS 

This work was supported by NASA Lewis Re- 
search Center through contract NAS 3-24105, project 
monitor B. Anderson. The author is indebted to 
B. McBride and S.T. Yu for their permission to use 
the equilibrium air package, and to Dr. J. S. Shuen 
for his helpful comments and suggestions. 

REFERENCES 

(l.j Pulliam, T. H., and Steger, J. L. "Implicit Fi- 
nite - Difference Simulations of Three Dimen- 
sional Compressible Flow.” AIAA Journal, Vol. 
18, No. 2, Feb. 1980, pp. 159- 167. 

[2.] Cooper, G. K. "The PARC Code: Theory and 
Usage.” AEDC-TR-87-24, October, 1987. 

[3.] Beam, R. and Warming, R. F. "An Implicit Finite 
- Difference Algorithm for Hyperbolic Systems in 
Conservation - Law Form.” Journal of Computa- 
tional Physics, Vol. 22, No. 1, Sept. 1976, pp. 
87-110. 


(4.) Pulliam, T. H. "Euler and Thin Layer Navier- 
Stokes Codes: ARC2D, ARC3D.” Notes for Com- 
putational Fluid Dynamics User’s Workshop, The 
University of Tennessee Space Institute, Tulla- 
homa, Tennessee, (UTSI Publication E02-4005- 
023-84), March 12-16, 1984, pp. 15.1-15.85. 

[5.] Jameson, A., Schmidt, W., and Turkel, E. "Nu- 
merical Solutions of the Euler Equations by Fi- 
nite Volume Methods Using Runge-Kutta Time- 
Stepping Schemes." AIAA Paper No. 81-1259, 
AIAA 14th Fluid and Plasma Dynamics Confer- 
ence, Palo Alto, California, 1981. 

[6.] Yu, S. T., McBride, B. J., Hsieh, K-C, and Shuen, 
J-S "Numerical Simulation of Hypersonic Inlet 
Flows with Equilibrium or Finite Rate Chem- 
istry." AIAA-88-0273, January, 1988. 

[7.J McBride, B. J., and Gordon, S. "FORTRAN 
IV Program for Calculations of Thermodynamic 
Data,” NASA TN D-4097, 1967. 

[8.] Moore, C., "Selected Tables of Atomic Spectra.” 
NSRDS-NBS 3, 1975. 

(9.) Moore, C., "Atomic Energy Levels.” Vol 1, 
NSRDS-NBS 35, 1971. 

(10.) Glushko, V. P., Gurvich, L. V., Bergman, G. A., 
Veyts, I. V., Medvedev, V. A., Khachkuruzov, G. 
A., and Yungman, V. S. "Thermodynamic Prop- 
erties of Individual Substances.” Academy of Sci- 
ences, Moscow, USSR, Vol 1, Parts 1 and 2, 1978. 

[11.] Balakrishnan, A., Davy, W. C. and Lombard, 
C. K. "Real-Gas Flows about Three-Dimensional 
Configurations.”, AIAA 83-0581, Jan, 1983. 

[12] Shuen, J. S., and Yoon, S. "Numerical Study 
of Chemically Reacting Flows Using an LU 
Scheme.” AIAA-88-0436, Jan, 1988. 

[13.] Gordon, S., and McBride, B. J. "Computer Pro- 
gram for the Calculation of Complex Equilibrium 
Compositions, Rocket Performance, Incident and 
Reflected Shocks, and Chapman-Jouguet Detona- 
. tions.” NASA SP-273 Interim Revision, 1976. 

[14.] Gnos, A.V., Watson, E. C., Seebaugh, W. R., 
Sanator, R. J. and DeCarlo, J. P. "Investigation 
of Flow Fields Within Large-Scale Hypersonic In- 
let Models.” NASA Tn D-7150, April, 1973. 


4 








P2 INLET: MACH 10 RAMP PRESSURES 
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Figure 8. Wall NO Mass Fractions, Moo = 10 Ramp Figure 11. Cowl Pressure Comparisons, P2 Inlet Ge 
Case ometry, Moo = 10, H = 100 Kft 


P2 INLET: MACH 10 EXIT TEMPERATURES 
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